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I. INTRODUCTION 



Bose-Einstein condensates (BECs) in trapped dilute gases exhibit interplay between quan- 
tum and nonlinear phenomena, since at zero temperature the mean-field Gross-Pitaevskii 
(GP) equation [l| for the order parameter applies with a good_accuracy. The macroscopic 

_and explained theoreti- 
One of the manifes- 



quantum coherence of BEC was demonstrated experimentally 
cally 4] with the use of the GP equation (see also the review Ref. 
tations of nonlinear dynamics in BEC is appearance of solitons, i.e. the self-localized "waves 
of matter", in quasi one-dimensional (cigar-shaped) condensates. 

Nonlinear phenomena in BECs, and the solitons in particular, are similar to that in 
nonlinear optics, for instance, in optical fibers j^. In both fields the nonlinear Schrodinger 
(NLS) equation appears at some level of approximation. Similar to optics, where bright 
and dark solitons are supported respectively by the focusing and defocusing nonlinearities. 



itons are now 



irl. while the 



in BECs the s-wave scattering length is the determining factor. Dark_so 
routinely observed in the condensates with repulsive interactions 0, IT 
condensates with attractive interactions allow for stable bright solitons 

The similarity with optics, however, does not go so far as the control over nonlinearity 
in the governing equations. Atomic interactions in BEC can be tuned at one's will by 
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application of magnetic field near the Feshbach resonance This opens a possibility 

to observe new phenomena in BEC when the atoms interact differently in different spatial 
locations of the trap due to the presence of an external field. Though such a setup has not 
been yet realized experimentally, it cannot be rendered as impossible at the theoretical level. 
For BECs constrained to lower dimensions, the Feshbach resonance proves to be sharp due 
to many-body effects, for instance, it is so in the two-dimensional condensate Similar 
result can be expected for the effectively one-dimensional BECs. The feasibility of control 
over the scattering length in BEC by the optical means is also being discussed 15| (see also 
Refs. ^Q). 

Control over the scattering length in a part of the condensate can be realized in the double- 
well trap with far-separated wells, with the gradient of the external field being localized in the 
region of the barrier where the order parameter is exponentially small. For such setup the 
well-known two-mode approximation applies, which considerably simplifies the governing 
equations and allows for the analytical study. The asymmetric double-well potential is 
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created by focusing an off-resonant intense laser beam near tlie center of a parabolic magnetic 
trap. Such setup is routinely realized in experiments (see, for instance, Ref. ^). 

The goal of the present paper is to understand the nature of the unusual bright solitons 
in the condensate trapped in an asymmetric double-well potential, when the larger share of 
atoms is localized in the well where the atomic interactions are vanishing or weakly repulsive. 
Such counterintuitive solitons, due to the usual association of bright solitons with attractive 
interactions, were recently numerically discovered 0]. However, the analytical form of the 
solitons and the domain of their existence were not addressed in the previous publication. 
In the present paper we analytically study the solitons and the domain of their existence. 

Nonlinear phenomena in BEC in the double-well trap have already attracted a great deal 
of attention. The coherent atomic tunneling between the two wells and the macroscopic 
quantum self-trapping phenomenon were theoretically predicted in Refs. jiol 



2G 



2l|. In 



Ref. 2l|, the predictions of the mean field approach based on the GP equation were also 
compared with the results from the full quantum model. The conclusion is that the char- 
acteristic time scale on which the quantum collapse and revival dynamics starts to play a 
significant role is large and increases with the number of BEC atoms, while on the smaller 
time scale the GP equation gives a quantitatively accurate description. Similar conclusions 
were derived recently in Ref. |22]. The general coupled-mode theory for the double- well 



trap with strong tunneling was developed in Ref. 



23| using the concepts of the nonlinear 



guided-wave optics. Nonlinear modes of the condensate in the double-well trap were used as 
the basis functions. Finally, the possibility of observing a chaotic dynamics of the conden- 
sate trapped in a double-well potential was the subject of Ref. 3]- An effective amplitude 
equation was derived, which exhibited a behavior resembling that of the Lorentz system. 

The paper is organized as follows. In the next section, complemented with the appendix, 
we give a self-consistent derivation of the coupled mode system for the asymmetric double- 
well trap with weak tunneling through the barrier. The general properties of the solitons 
in the coupled-mode system are summarized in section IIII Al In section IIII B| we derive 
an approximate analytical form of the soliton solutions in the case of vanishing interactions 
in one of the condensates. The exact boundary of the soliton bifurcation form zero, which 
defines the domain of existence of the stable solitons with the larger share of atoms in the 
repulsive condensate, is found in section lTlI CI Section|TVlcontains discussion of the predicted 
phenomena. Finally, all numerical simulations in the present paper were performed with the 
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use of the spectral collocation methods 
be found in Ref. Q- 



25| known for their high-accuracy. More details can 



II. REDUCTION OF THE GP EQUATION TO THE COUPLED-MODE SYSTEM 

Consider BEC trapped in an asymmetric double-well potential and assume that each of 
the wells is weakly confining in the x-direction and strongly confining in the transverse plane 
rj_ = (y, z). Let the transverse trap be parabolic in the z-direction and have the double- well 
shape in the ?/-direction, see fig. [H 

The zero-point energy difference between the two wells (shown by the dashed lines in 
fig. ^ can be managed, moving the position of the barrier-generating laser beam. We can 
neglect the small variation of the trap potential in the x-direction. The corresponding 
Gross-Pitaevskii (GP) equation for the order parameter \Ef of the condensate then reads 

«^ = l~(l, + Vi ) + V„.(rx) + *, (1) 



dt \ 2m ydx"^ 
with the external potential 

V„. = ^fe= + 7V) + V.exp{-(|^}. Vs>0. (2) 

Here the parameter 7 accounts for the trap asymmetry in the transverse plane, while the 
parameters yo and a give the position and the width of the barrier, respectively. 

The double-well potential Vext is assumed to have two quasi-degenerate energy levels: 

El — Eq <^ E2 — El, (3) 

fig. H can serve as an illustration. There is a natural choice of the basis wave functions in 
the degenerate subspace, say ipu and ipv, where each one is localized in just one of the wells. 
Such basis is given by a rotation of the wave functions for the ground and the first excited 
states: 



The choice of the parameter x is based on the fact that the quotients of the absolute values 
of the eigenfunctions V'o(ri) and V'i(^±) at the extremal points ?/„ and ?/_,_ (in the two wells) 
are in approximate inverse proportionality relation 

^i(y-) _ ^o(y+) 



^0(1/-) ^\{y+) 



(5) 
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FIG. 1: Illustration of the double-well: the parabolic trap potential in the y-direction modified by 
a Gaussian barrier. The (horizontal) solid lines indicate the numerically computed eigenvalues in 
the combined potential. The dashed lines show the bottom energies in the two wells of the trap. 
We use dimensionless units. 

(the positions of the extremals of the wave functions slightly deviate from each other and 
from the minima of the trap; these deviations we neglect). Equation holds for the basis 
functions of the degenerate subspace for a general asymmetric double-well trap with weak 
tunneling through the barrier (see the appendix). For the double-well of fig. ^the quotients 
of the eigenfunctions shown in fig. El have the values 3.1 at the left and —2.9 at the right 
minima of the trap. We set h = il)i{yJ)/il)Q{yJ). For this choice of x the wave functions 
'ipuiy) and Tpviy) defined by equation (jH) are localized in the left and right wells, respectively, 
see fig. 121 

For not too large number of atoms (see condition (jiup below) one can neglect the higher 
excited levels of the trap and approximate the solution to the GP equation (P) as follows 



This approximation allows for reduction of the GP equation to a system of linearly coupled 
equations (equations (|7|) below) for the projections f^^v of the order parameter on the local- 
ized basis, due to approximate decoupling of the basis functions in the nonlinear term. For 
the double-well given in fig. Q for instance, we obtain the following nonlinear cross-products: 



(6) 
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FIG. 2: The ground and first excited state wave functions for the double-well trap in fig. ^ (dashed 
lines) and the wave functions localized in the complementary wells (solid lines) constructed as 
indicated in the text. All units are dimensionless. 

Substituting expansion (jHl) into equation ((T)), multiplying by the localized wave functions, 
integrating over the transverse variables y and z, and throwing away small nonlinear cross- 
terms we arrive at the coupled-mode system of Refs. 
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but with the kinetic terms: 



ih 



ih, 



dfu 
dt 

dfv 



h' d'fu 

2m dx^ 



(7a) 



dt 2m dx^ 

The zero-point energies Eu and By and the tunneling coefficient K are defined as follows 



(7b) 



(8) 



(9) 



here the z-dependence of the wave functions ipu and ip^ is given by the ground state wave 
function (the Gaussian) in the 2;-direction (note that by changing the sign of either ip^ or 
V'w, if necessary, one can always set K > {]). Below, we assume that the atomic interaction 
in the f-condensate is externally modified. 

In the above analysis we have neglected the effect of nonlinearity on the dynamics in 
the transverse (y, 2;)-plane assuming a small contribution from the atomic interaction as 



compared to the energy gap E2 — Ei of the trap. This leads to the following condition for 
each of the two parts of the condensate: 

^ « 1, (10) 

where M is the number of atoms, is the scattering length, and I is the length of the 
condensate. Equation ()10|1 is the usual condition of applicability of the one-dimensional 
NLS equation for the cigar-shaped condensate. The condition of applicability of the classical 
mean-field approach is <^ 1 and A/" ^ 1 (see, for instance, Ref which is satisfied 

for large number of atoms and condition ()1U|) . Finally, the general applicability condition 
for our two- mode approximation reads K <^ hio^ (see the appendix). 

In the calculations it is convenient to use a dimensionless form of system (|7j) with the 
following dimensionless variables: 

c ^ . ( ^ 2{E,-E^) 2K a, 

\muj±J 2 nuj± nuj±^ \au\ 

Here au {clu < 0) and are defined as follows: au,v = d\ J asip^,v^'^^-L with being the 
s-wave scattering length (different in the two condensates). Introducing the dimensionless 
order parameters for the condensates in the two wells, 

«(e,r) = e^^(87r|a„|)V2/., ^^^^^^ ^ ^,^^g^|^^|^i/2^^^ ^^2) 

we arrive at the coupled-mode system suitable for the numerical study: 

du d'^u . ,2 / s 

+ + \u\ u + Kv = (13a) 



i— + - {e + a\v\^)v + Ku = (13b) 



. dv d'^v 

The number of particles in the coup led- mode system (fT^. defined a.?, = /d^|Mp and 
A^^ = / d^|f p, is related to the actual number of atoms N as follows: Nu,v = {.d^/STi\au\)Nu^^ 
(the coefficient d±/'&Ti\au\ is of order 10^ — 10'^). Note that the quotient of the number of 
atoms in the two condensates does not change under the scaling transformation. Below we 
will use only the quantity Nu^v referring to it as "the number of atoms" for simplicity. 

Finally, we note an important difference between our approach and that of Ref. [2^. 
In the latter work the nonlinear modes are used as the basis functions, what results in an 
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effective coupled-mode system with nonlinear cross-terms, whereas our system ()13p consists 
of linearly coupled equations. Our coupled-mode system, however, applies only in the case of 
weak tunneling through the barrier. But precisely in such setup the experimental realization 
of control over the scattering length in a part of the condensate may be feasible. 



III. SOLITONS IN THE TUNNEL-COUPLED CONDENSATES 

A. General properties of the soliton solutions to the coupled-mode system 

We are interested in the stable stationary soliton solutions to the coupled-mode system 
(jl3|) mainly for a > 0, i.e. for a repulsive-attractive pair of condensates. Setting u = 
e^^'^^U{C,) and v = e~^'^^V{^), where fi is the dimensionless chemical potential and the 
functions f/(0 ^-^e real, we obtain the corresponding stationary system for the 

soliton profiles: 

I^U + % + U^ + kV = 0, (14a) 



{fi-e)V + V^^-aV^ + kU = 0. (14b) 

The two-mode approximation requires that k <^ 1 and \e\ <^ 1. On the other hand, system 
()14|) admits the following scaling invariance: 

U = \/kU, V = \fKV , ^ = —;=, a = Ku, e = ne, (15) 



where the new variables with tilde satisfy the same system with k = 1. For this reason we 
will keep e and k arbitrary. The scale-invariance ()15|) can be used to verify the analytical 
solutions. 

We are interested in the so-called fundamental solitons, i.e. when the functions U{^) and 
V{C,) are nodeless and decreasing exponentially as |^| oo. The exponential decay rate k, 
U ~ e~'^^ and V ~ e~^^ as oo, is given by the dispersion law of the linearized system. 
One obtains the following two branches: 



£2 „ e 



ki,2 = - /ii = - - y J + ^2 = - + y J + (16) 

Note that for arbitrary k, and e the boundary chemical potentials satisfy the inequality 
fii < < fi2- They coincide with the two quasi-degenerate energy levels of the coupled- 
mode system when the nonlinearities are negligible. 
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For a > it can be shown that the two dispersion laws correspond to two possible classes 
of the fundamental soliton solutions: UV > or UV < 0, below called the in-phase and out- 
phase solitons respectively. Indeed, the amplitudes Uq = U{0) and Vq = V{0) are bounded 
as follows [isj ]: 

< -2/i, if UV > 0; aV^ < 2(/x - e), if UV < 0. (17) 

Noticing that the conditions kf > and > demand respectively fi < min(0, e) and 
fi > max(0,£:) we see that the in-phase solitons correspond to the fci-branch, while the 
out-phase ones pertain to the /c2-branch. Numerical calculations confirm this prediction 

0. 

We will discard the out-phase soliton solutions from consideration, since they bifurcate 
from zero at a higher chemical potential /i2 (fTBj) . Therefore, it is clear that they can be stable 
only in e xcep tional cases (they do not satisfy stability criterion (j?H) given below, consult 
also Ref. [22|). Below by the soliton solution we will mean the in-phase soliton. 

When the condensates are attractive, i.e. for a < 0, there is a class of the sech-type 
soliton solutions which generalizes the well-known symmetric solitons in the standard 



model of the dual-core optical fiber 
soliton solution is given as 



see also the review Ref. 



). The sech-type 





sech (^^^ , An = (-/i//t - (-a)-'/') . (18) 

It exists in the special case, when e, k, a, and n are related as follows: 

£ = K [(-a)^/2 - i-a)-^/^] = h{K, a), /i < /ii = -K{-a)-^/^. (19) 

Evidently, the exact sech-type soliton solution does not exist for a > 0. 

There is another branch of the in-phase solitons, not known analytically, which bifurcates 
from the sech-type soliton branch (fTH|). in the domain of their existence, at the point 

m = t r-^ , ^ = Kf^, «)• (20) 

6y/—a 

For other parameter values the sech-type solitons (jl8|) and the bifurcating in-phase solitons 
coexist on separate curves = N{fi) (here is the total number of atoms in the 
condensate). At least one of the two soliton branches must have continuation for positive 



values of a, since we have found such in-phase sohtons numerically in Ref. (we have 
found exactly one branch for a > 0). 

n 

It is important to know the stability properties of the soliton solutions. In Ref. |18|| we 
have shown that for a > the Vakhitov-Kolokolov criterion for the soliton stability 

dN 



9/i 



< 



(21) 



can be applied (consult also Refs. 30,|31|). The stability property of the in-phase solitons 
for a < is affected by bifurcations (i.e. criterion (j21|) does not apply for all /i when a < 0). 
Details can be found in Ref. jisl . 



B. Analytical soliton solutions for a = 

For a = the self-bound states, solitons, are possible due to the balance between the 
cubic nonlinearity in the attractive condensate and the linear dispersion in both condensates. 
Thus we can use the following strategy to find an approximate expression for the solitons. 
Equation ()14b|) . which can be cast as 



AV=(^1- /3-^^ V = K/3U - a(3V\ /3 = - yu)-\ 



(22) 



can be explicitly resolved with respect to V for a = using the method of the Green 
function. It was already mentioned that for the in-phase solitons fi < e and hence (3 > 0. 
Thus the operator A is positive definite and invertible and the unique Green function reads 
G{^) = exp{-|^|/V^}/(2v^). Therefore the solution to equation ^ for a = is 



oo 



V{0 = ^ / e-^-^Uir])dr]. (23) 



-oo 



2 

Solution (j^Hj) can be expanded in the infinite series with respect to the derivatives of U{C,): 

ViO = (uiO + f;rf/('")(0 ) , (24) 



n=l 



which can be most easily verified by the formal inversion of the operator A. The crucial step 
in the derivation of the approximate soliton solution consists in retaining only two terms 
from the series ()24j) when it is substituted into the equation for U{C,), equation ()14a|) . For 
such a reduction the following condition must be satisfied 

> fi;2/3«+if/{2n)^ n>2. (25) 
10 



Thus, under condition the coupled-mode system (fT^ reduces to the canonical NLS 

{H + K^(3)U + {l + K'^(3'^)Ui:i: + U^ = 0. (26) 
In this approximation the w-component of the solution is given by the sech-type soliton: 



U = Asech(fcO, A 



2{ni - n){n2 - /i) 
e — fi 



1/2 



k 



(£-/i)(/il -^)(^2 -/i) 

{e — /i)^ + 



1/2 



• (27) 



Using the identities /ii/i2 = —i^"^ and e = yUi + /i2 it is easy to check that the amplitude 
A satisfies the first inequality in equation ()17p for all values of fi. However, the analytical 
decay rate k is not the same as ki in equation (fT^ . This is the price we pay for adopting 
the approximation. In fact, as either /i ^ /ii or ^ (recall that ^ < fii < 0) the 
decay rate k approaches ki. 

k = ki[l + C(/ii -//)], k = ki[l + 0(/xi//x)] . (28) 

The applicability condition for the analytical approximation is then satisfied at least in these 
two limits. Indeed, condition demands that A^ ^ n'^f^^k'^A, where we have estimated 
~ kU, and after a simple transformation we obtain 



2 _ (/ii -/i)(/i2 -^) , {e-fif . . 



Moreover, comparison with the numerical solution shows that the analytical approximation 
captures all essential features of the number of atoms as function of the chemical potential 
for all values of fi (see fig. El below). 

The f -component of the soliton, given by equation ()23|) . cannot be found in the explicit 
form. Nevertheless, all the necessary physical quantities, for instance, the number of atoms 
Ny{fi) and the energy of the soliton, can be found in closed form. The shares Nu and Ny of 
the BEC atoms trapped in the two wells read: 

iV. = 4((^^-^^(^^-^^^3-^^' + "']y'', (30) 

2 (/^i -/^)(/^2 -m) (2n + l)a + 3 

where a is defined in equation (jSHl)- 

The deformation of the curves N^^fi) and Ny{fi) given by equations (jHUI) and (jSD) with 
variation of the zero-point energy difference e is illustrated in fig. El (we remind that the scaled 

11 



number of atoms = (87r|a„|/(i_L)A^ is used). The qualitative behavior of the analytical 
curves is in excellent agreement with the numerical results of Ref. (the local maxima of 
the analytical curves Nu^^{ii) are narrower and 10% higher than the maxima of the numerical 
curves) . 
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FIG. 3: Behavior of the number of atoms in the two wells of the trap vs. the chemical potential 
with variation of the zero-point energy difference between the wells (all quantities pertain to the 
dimensionless system Here the dotted lines correspond to N^^ the dashed lines to A''„, and 

the solid ones to the total number of atoms = Ny^ + Ny In panel (d) only the part of the picture 
about the local maximum is shown (the full picture is similar to that in panel (c)). 



Note the striking sensitivity of the total number of atoms and that of the fraction N^/Nu 
at the local maximum with respect to variation of the zero-point energy difference e. From 
figs, infc) and (d) it is seen that a small variation of e (with e being on the order of a small 
tunneling coefficient k) results in a large variation of both the total number of atoms and 
the relative share of atoms in the w-condensate. 
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The solitons corresponding to the two stars in fig. Efd) are shown in fig. lU Sohton-(l) is 
unstable while soliton-(2) is stable. The latter one represents the unusual self-bound state 
of the system, where almost all atoms are gathered in the noninteracting condensate. 
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FIG. 4: Unstable (1) and stable (2) soliton solutions to the dimensionless system (jl3|) with the 
same total number of atoms which correspond to the two stars in panel (d) of fig. IHI 



From fig. El it is seen that the soliton solution bifurcates from the zero solution at /i = /^i. 
In fact, for a = 0, the total number of atoms is indeed always equal to zero on the boundary of 
the soliton existence. For ^ ^ fii the sum on the r.h.s. of equation (jHlll can be approximated 
as follows 

^2{fii-fi)-'^' + 0{{fi,-f,Y/'}, (32) 



(2n + l)a + 3 
^ [(2n + l)a + l]3 



n=0 



where have used the Euler-Maclaurin summation formula and the definition of a ()29j) . Hence, 
the number of atoms is proportional to {fii—jJiY^'^ as /i — /ii and the bifurcation is continuous 
for a = (in contrast, the number of atoms has a finite limit a.s fi —>■ fii for a > and e = 
see the next section). 

From equations (jHOJ-lISS) one may deduce the ratio of the slopes of A^„(yu) and Ny{fi) as 

AT I,, J , /TTJoWn^ - r/9 

(33) 



hm — 



^(e/2)2 + k2 - e/2 



/i2 ^(e/2y + K^ + e/2' 
This is the quantitatively accurate result, since in this limit condition ()29p is well satisfied. 
Formula ()33|) predicts, for instance, that almost all atoms are gathered in the non-interacting 
condensate for |e| k and e < 0. Thus, to observe such unusual solitons the zero-point 
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energy difference e should be decreased sufficiently below zero. Physically, it is clear that 
the energy decrease due to the attractive interaction in the w-condensate with increase of 
should be compensated by a lower zero-point energy in the w-well, if one wants to increase 
the ratio N^/Nu- 

Condition (j29|) breaks down about the local maximum of the curve N = N{fi) (see 
fig-Etd)) and one has to use the whole infinite series expansion f)24|) for quantitatively correct 
description of the soliton in the system. Therefore, the condensates are strongly coupled 
when almost all atoms are trapped in the f-well. We will see below that the unusual solitons 
also appear for < a <^ 1 and Scr < £ < 0, with almost all atoms being contained in the 
repulsive condensate strongly coupled by the quantum tunneling to the attractive one. 

For ^ \fii\ we obtain the following estimates for the number of atoms in the two wells: 



N,.=AV^+0{\^^\-'/% ^ 



0.712^ + O 



(34) 



where we have used a = 1 to evaluate the infinite sum in equation (jM]) (since a = 1 + 
(9(yU^^)). First of all, almost all atoms are gathered in the attractive w-condensate. Second, 
in this limit the w-component of the soliton solution asymptotically satisfies equation ()14a|l 
for K = 0, since the amplitude and decay rate of soliton (j77|) asymptotically approach 
those for the soliton solution of a single NLS equation: A = \J—2[i + (9(|/i|~^/^) and k = 
+ (9(|/i|~^/^). Thus, for large negative jj the coupled-mode system predicts weak 
coupling between the two condensates. 

To determine the ground state of the system let us find the energy of the soliton solutions. 
The GP equation corresponds to the energy functional 



dV<!— |v^l^ + Vexti^r + |l^l 



(35) 



In the coupled-mode approximation we have the following expression for the energy 

£ = EM+^^M, (36) 
where M is the actual total number of atoms and H is the Hamiltonian for system (fTT^ : 



H 



du 



dv 



I I 1 2 / * I *\ 

+ e\v\ — K[uv + vu ) 



\u\ 



+ 



a\v\ 



The Hamiltonian H evaluated on the soliton solution reads 



sol 



(/il - /x)^(/i2 - /i)' 

{e — /i)2 + 



1/2 



{e - lifl^ 



4k- 



. (/ii - fi) {fi2 - fJ') (2£ - 3/i) 
{e - /i)V2 



(37) 



(38) 
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where S{fi) = X^^^o [(2n+i) "+ip ^ ~ if^l^I- Taking into account that 2e > max(3yU, 5yu), 
since /i < /ii, it is easy to see that ifgoi is always negative. In the coupled-mode approxi- 
mation, expression ()H8j) is quantitatively correct for /i ^ /ii and \^\ ^ and for other 
values of fx captures the qualitative behavior of the energy. Note that H scales as k^/^. 
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FIG. 5: The energy vs. the total number of atoms (for the dimensionless system ()13() ') for two values 
of the zero-point energy difference: £ = (the dashed line) and e = —k (the sohd hne). The sohd 
hne has a very narrow twist without self-intersection, which occupies the interval 8 < N/ ^/k < 12. 

Lowering e below zero forces the energy to develop a twist, see fig. El The twisted part 
of the curve widens by decreasing e, while its left endpoint moves towards = 0. (The 
curve H = H{N) has the shape similar to that of /i = fi{N) after the following substitution 
is applied to the latter: (/ii — fi) —H.) It follows that, for almost all values of e and the 
total number of atoms, the energy is minimized on the solutions with large negative /i (in 
fig. El the latter solitons correspond to the lower part of the solid line below the twist, while 
the unusual solitons pertain to the twisted part). However, according to criterion (PT|) the 
solitons with the larger share of atoms in the f-condensate are stable (precisely, the solitons 
which pertain to the upper part of the twist in fig. thus they can be experimentally 
observed. 

Finally, let us check condition (fTIUl . Estimating the condensate length as i/d± ~ k, we get 
that in the dimensionless variables the one- dimensionality condition for the -u-condensate is 
kNu <^ 1. Using the explicit expressions (P7|) and (^01) we arrive at (/xi — /x)(/i2 — f^) ^ s — ii 
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which is satisfied for /i close to /xi, i.e. in the domain of the unusual solitons. Similar, for the 
f -condensate we have fc|a|A^^ <^ 1, which is satisfied for small a in the same range of /i. On 
the other hand, the soliton solutions with large negative /i may violate the one-dimensionality 
condition (it is violated in the limit ^ It is a refiection of the fact that for large 

enough number of atoms the transverse degrees of freedom start to play a significant role 
for the weakly coupled condensates, i.e. for the ground state in the system, leading to the 
collapse instability at some kN = f{k) (for details on collapse consult Ref. js^). However, 
for the number of atoms below the collapse instability threshold the ground-state soliton 
of the coupled- mode system is stable (but may have a modified shape), while the unusual 
solitons are effectively one-dimensional, since they always satisfy condition (fTUj). 

C. Solitons for a > and the boundary of the bifurcation from zero 

We have seen that for a = the soliton solutions bifurcate from the zero solution at 
fi = fii and the unusual stable solitons appear close to the bifurcation boundary. For a > 
the competition between the zero-point energy difference e and repulsion in the f-condensate 
creates a lower boundary = &(^, d) for the bifurcation from zero. Thus for a > the 
unusual solitons appear only for e > Existence of such a boundary was discovered 
numerically in Ref. |18| for a = 1. 

It turns out that ecr{k,a) can be found analytically in the exact form. To this goal we 
suppose that the solitons do bifurcate from zero and find the value of e for which we have a 
contradiction. Setting e = — fi, we assume that U = 0{eP) and V = 0{e'^) as e 0, for 
some p,q > 0. Using the estimate d^/d,^^ ~ A;^ = e and P = (/i2 + e)^^ = 0{1) in equation 
()22|1 we immediately obtain q = p. Note that the terms containing U in equation p4aj) have 
the following powers of e: /if/ = C(e^), U^^ = C(e^+^), and = 0{e^^), where we have used 
the above estimate for the second derivative. Therefore, irrespectively of the value of p, it is 
clear that if an expression for V contains all the terms with powers of e up to max(p + 1, 3p), 
it then contains all the necessary terms we need to derive the soliton solution in the lowest 
order of e. 

The most general expression for V that explicitly contains all the powers of e up to 
max(p + 1, 3p) is as follows 

V = F{U) + CU^^ + 0{e'), s = imn{3p + l,p + 2), (39) 
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where the function F{U) and the constant C (of order 1) are determined by solution of 
equation (j^^ in the required order of e. We will need the expression for V up to the terms 
of order e^^ only. We obtain 



V = K(3U - aK-'p^U' + 3a'K''(3^U'' + Kl3'U^^ + 0{e'), r = min(3p + l,p + 2, 7p). (40) 

Here the coefficient P = {^2 + also should be expanded into powers of e. 

Substitution of expression pUjl into equation ()14aj) leads to the following equation for the 
M-component of the soliton: 



-e 



l-^]u+(l-^)u,,+ (l-'4]u' + '-^U^ = 0{en. (41) 



/i2/ V /^2/ V / /^2 

From equation (j4ip it is quite clear that the bifurcation from zero is forbidden if the 
coefficient at the cubic nonlinear term is negative (we remind that ^2 > /Wi), i-e. if the 
lowest-order nonlinearity is repulsive. In the latter case equation ()41|) does not contain 
bright soliton solutions in the limit e 0. Hence, the exact boundary of the bifurcation 
from zero is given as /i2(£cr, f^) = ko^/^. Using the definition of /i2 ()16|1 . we obtain the critical 
value of the zero-point energy difference 

e^^ = n (a^/^ - a'^'^) . (42) 

There are no in-phase solitons bifurcating from the zero solution below the boundary 

The properties of solitons close to the bifurcation boundary /ii are different in the two 
cases e > and e = Scr- In the former case equation ()41|) reduces to the canonical NLS 
equation (in this case p = 1/2) and the corresponding soliton solution reads 



2nl{n2 -/ii) 

fj,2 — O'f^'^ 



nl/2 



(43) 



The u-component of the soliton, obtained from equation ()40|) . in the same order in e reads 
V = -U + 0{€'/^) = ^v^sech( + (^(e'/'). (44) 

/^2 /^2 

It is straightforward to verify by integration that the number of atoms corresponding to the 
soliton solution pHjl - ljUjl goes to zero as e ^ 0: 

iV„ = ^4^^^v^+0(e^/^), N^ = ^-^4^^^^^e + 0{e^l-). (45) 

/i2 — OK^ /i2 — 0./^ 
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We conclude this case noticing that for e > the sohtons bifurcate from the zero solution 
continuously. 

It is interesting to note that if a < and e = h{K,, a), i.e. as in equation (|19p. the soliton 
solution given by equations (j43|) - (j44j) coincides with the sech-type soliton solution (jl8p (in 
this case fxi = — k(— a)~^/^ and jj,2 = k(— a)^/^). Therefore, in this special case, one can 
think of the family of solitons ()43|) - P^ as a continuation of the sech-type soliton solutions 
(UHl) for a > 0. 

In the critical case, e = the cubic term is absent from equation ()4H) . The effective 
equation for the w-component of the soliton is the quintic NLS equation (now p = 1/4). The 
corresponding soliton solution reads 



U = Boe'/' sechi/2(2v^O + (^(e'^'), Bo 



fi:(l + ai/2) 



,3/4 



1/4 



(46) 



Here we have used the relations /i2 = na^^'^ and /ii = —Ka which follow from equation 
(fTH|) for 6 = Ecr- The w-component of the soliton is given by 



V 



-U + 0{e 



3/4 N 



(47) 



Though soliton solution pH|) - (j77j) bifurcates from the zero solution, since both components 
tend to zero in the limit e — 0, the corresponding number of atoms now approaches a finite 
limiting value. Indeed, we obtain: 



IT 

2 



k{1 + a 



l/2\-| 1/2 



+ 0{e 



3/2 ^ 



IT 

2 



k(1 + a 



l/2^ 



1/2 



+ 0{e 



3/2N 



(48) 



The critical bifurcation can thus be rightfully called a discontinuous one, since the total 

number of atoms must be greater than some threshold value for the system to have a soliton 

solution of vanishing amplitude. 

In fig.iniwe show the numerically computed behavior of the number of atoms vs. chemical 

potential for a = 0.001 in the two cases: e > (panel ( ) and on the boundary 

of the bifurcation from zero e = (panel (b), = —5.45k). 

It is seen that the curves shown in fig. El^b) indeed cease at the finite values Nu = 21 and 
= 670 (according to equation (j48j) ) on the boundary fi = fii of the soliton existence. The 

corresponding soliton solutions are unstable as /i — > /ii, as they correspond to the positive 

slope of the curve N = N{fi) for the total number of atoms. Solitons corresponding to 

fig- Efa), on the contrary, are stable in this limit. 
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FIG. 6: Bifurcation of the soliton solutions for two values of the zero-point energy difference: 
e = —5.30k (panel (a)) and ^cr = —5.45k (panel (b)). The dimensionless units of system ((T^ are 
used. Here a = 0.001. In both panels: the dotted lines correspond to N^, the dashed to A',;, and 
the solid ones to the total number of atoms. 

Moreover, for e < Ecr, for values of the chemical potential close to the boundary fii 
there are soliton solutions which represent the unusual self-bound states of BEC, where 
almost all atoms are gathered in the repulsive f-condensate. The unusual solitons have the 
shape similar to that of soliton-(2) of fig. EJand appear for small values of the parameter a. 
Indeed, their appearance is conditioned by e < at least, according to formula (jHSj), while 
the condition a < 1 is necessary for a negative lower boundary Obviously, the stability of 
such solitons is derived from the strong tunnel-coupling to the attractive u-condensate and 
the atoms can be thought as existing in the "superposition" of the repulsive and attractive 
pairs. 



IV. CONCLUSION 



We have studied the solitons in BEC in a double-well trap with far-separated wells, when 
the atomic interaction is externally switched from attractive to repulsive in the condensate 
trapped in one of the wells. 

First of all, there are soliton solutions with almost all atoms contained in the attractive 
condensate, which correspond to weak tunnel-coupling between the condensates in the two 
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wells of the trap. In this case the attractive condensate serves as a spatially inhomogeneous 
source of atoms for the repulsive one, thus tailoring the repulsive condensate to the sohton 
form. 

However, there are unusual stable solitons in the system, when almost all atoms are 
contained in the repulsive condensate. It is found that such soliton solutions appear for 
weak repulsion when the zero-point energy of the well with the repulsive condensate is lower 
than that of the attractive one but is above some critical value depending on the tunneling 
and interaction coefficients. The unusual solitons describe the condensates strongly coupled 
via the quantum tunneling, notwithstanding the fact that the tunneling coefficient is small. 
This apparent paradox is resolved by noticing that the tunneling probability is proportional 
to the number of atoms and that the total number of atoms increases with increase of the 
relative share of atoms in the repulsive condensate. 

The ground state of the system corresponds to the soliton solution with the larger share 
of atoms contained in the attractive condensate. However, the solitons with the larger share 
of atoms in the repulsive condensate are stable, hence they can be experimentally observed. 
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APPENDIX: WAVE FUNCTIONS WITH THE LOCALIZATION PROPERTY 

Here we show that the inverse-proportionality relation (jSj) holds for the eigenfunctions 
corresponding to the quasi-degenerate energy levels in the double-well trap with weak tun- 
neling through the barrier. 

The weakly asymmetric double- well potential can be cast as V = Vs{y) + Va(?/), where 
'^si~y) = '^s{y) is the symmetric potential and Va(— = — Va(i/) is the perturbation. 
Under assumption (jH)) the eigenfunctions of the full Hamiltonian 

H = Hs + VAiv) ^ -^-f^ + Vs{y) + VAiv), (A.l) 
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corresponding to the shifted degenerate energy levels are given by a rotation of the wave 
functions |0) and |1) corresponding to the degenerate energy levels of the unperturbed Hamil- 
tonian Hs'- 

where AE" = Ei — Eq and A£^ = Si — So is the energy spitting in the unperturbed and 
perturbed trap, respectively, and K = (0|Va|1). The condition for equation ()A.2|) is K 
hxj±, where hw±_ serves as an estimate on the energy gap E2 — Ei in the double- well trap. 

Equation ()A.2|) guarantees that if the inverse proportionality relation © holds for the 
symmetric double-well it also holds for the asymmetric one. This follows from a simple 
calculation: 

fa+|o)+-^ _ {y+\^i) . . „x 

^M>_1 {y+\SS ^ ' 

where the supposed inverse proportionality relation (0) for |0) and |1) have been used. On the 
other hand, under the assumption that the tunneling coefficient is weak the eigenfunctions 
of the symmetric double well are well-known combinations of the eigenfunctions for the 
two far-separated wells: |0) = + |r))/v^, |1) = {\tj — \r))/\/2. Then the relation 
T^^TTT ~ — ~ 1 is obvious due to localization of the wave functions \t) and |r) in the 
respective wells and the symmetry. This concludes the proof of equation Q. 
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